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INTRODUCTION  AND  SUMMARY 

Asymptotic  approximations  of  optimal  control  laws  are  determined  for  a 

class  of  multivariate  dynamic  systems  whose  state  components  are  either 

slowly  varying  (i.e.,  parameters)  or  can  be  tracked  arbitrarily  well  by  making 
the  controller's  state  measurements  sufficiently  accurate.  This  control 
optimization  problem  would  have  the  standard  linear-quadratic-Gaussian  form 
except  for  small  terms  in  the  measurements,  that  are  bilinear  in  the  control 

and  parameter  variables,  and  with  respect  to  which  the  control  law 

approximations  are  asymptotic.  Also,  the  measurement  noise  is  small  in  a 
certain  relative  sense,  which  gives  this  control  problem  special  properties. 

The  bilinear  measurement  terms  give  rise  to  a  rapidly  fluctuating  term  in  the 

optimal  control  law — a  linear  function  of  the  state  of  a  multivariate  linear 

system  driven  by  a  Kalman-filter  innovation  variable. 

This  rapidly  fluctuating  control  component  is  essentially  zero-mean, 

even  conditioned  on  all  but  very  recent  data,  and  so  has  little  current  effect  on 
the  system  dynamics.  It  is  apparently  an  example  of  the  probing  phenomenon 
in  optimal  control  laws  identified  by  Feldbaum  (Reference  1),  which  is  current 
control  effort  expended  to  reduce  uncertainty  in  the  state  variable  (including 
the  parameters  here)  in  order  to  improve  system  performance  in  the  future. 
A  voluminous  amount  of  literature  exists  on  this  general  subject.  In  work 
related  to  the  particular  topic  here,  specialized  methods  have  been  u^ed  to 

obtain  sharper  results  for  a  particular  case  related  to  homing  missile  guidance 
(Reference  2).  Speyer  and  Hahn  have  derived  similar  asymptotic 
approximations  of  optimal  control  laws  for  systems  with  bilinear  terms  in  the 

dynamics  (Reference  3),  but  the  "parameters"  there  are  not  static  or  slowly 
varying  and  the  measurement  noise  is  not  small. 

The  results  here  are  obtained  from  an  approximate  dynamic 

programming  analysis  in  which  the  state  variable  includes  the  departures  of 
certain  Kalman-filter  covariance  quantities  from  their  local  averages,  and  the 

analysis  uses  the  fact  that  these  departures  are  small  because  of  the  small 
measurement  noise.  It  is  then  shown  how  these  results  apply  to  other  cases 
that  arise  from  analyses  of  noise-induced  departures  from  nominal  behavior 
in  a  more  general  class  of  optimal  control  problems.  Typically,  these  other 
cases  have  (relatively  small)  quadratic  terms  in  the  system  dynamics  and  state 
measurements,  and  cubic  terms  in  the  performance  criterion.  The  analysis 
here  is  not  at  a  mathematically  rigorous  level,  although  the  constructions 
developed  might  be  useful  in  a  more  ambitious  treatment  of  that  sort. 
Expressions  denoting  ordinary  differential  equations  with  white  noise  terms 
should  be  understood  as  the  formally  corresponding  stochastic  differential 
equations  in  the  Ito  sense  (Reference  4  and  5)  if  a  rigorous  interpretation  is 
desired. 
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Unless  otherwise  stated,  lower  case  letters  denote  (real  finite¬ 
dimensional)  column  vectors  and  scalars.  Matrices  are  denoted  by  capital 

Roman  letters.  ATdenotes  the  transpose  of  a  matrix  A,  and  lr(A)  its  trace  if  A  is 
square.  It  will  also  be  convenient  to  make  use  of  three-way  matrices,  which 
are  always  denoted  by  capital  Greek  letters  here.  For  continuity  of  notation, 
the  following  definitions  are  adopted  for  such  a  three-way  matrix  r,  with 
vector  x  and  matrices  A  and  B  of  compatible  dimensions  and  with  repeated 
indices  denoting  summation: 


(rx)ij  =  rijoxo 

(matrix) 

(Ax7)ijk  =  Aijxk 

(three-way 

matrix) 

(Ar)ijk  ~  -^a^ajk 

(three-way 

matrix) 

(FB)ijk  =  FjjjjBgk 

(three-way 

matrix) 

)ijk  ~  ^"jki  )ijk  ~  ^kji 

(three-way 

matrix) 

[trcrjij  =  roia 

(column  vector,  when 

applicable). 


T 

With  these  definitions,  the  expression  FArBDxx  is  fully  associative.  Many 
other  consequences  are  obvious.  Some  useful  but  less  obvious  properties  are 

tr(r'x)  =  [tr(r)]Tx 

A  tr(r')  =  tr[(Ar)'] 
tr(AT)  =  tr(r  A ) 

(rB)'  =  bV  and  (Ar)"  =  r"AT 
(ArB)T  =  btgtat 

(r’x)AT  =  (Ar)'x  and  (r"x)B  =  (rB)”x. 


PROBLEM  STATEMENT 

The  more  limited  case  actually  analyzed  here  is  that  of  a  multivariate 
system  with  dynamics  of  the  form 

|  x  =  Fx  +  Gu  +  w  (  1  ) 

(  Q  =  h  w2,  (  2 ) 
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a  controller  of  which  receives  a  vector  measurement  of  the  form 

z  -  Hx  +  h  H  tr(r"8uT)  +  v  (3) 

and  selects  the  control  vector  u  at  each  time  instant  t  >  0  (with  the  convention 
that  Equation  3  is  evaluated  slightly  before  it  is  used  to  change  u).  The  initial 
state  part'tions  x(0)  and  0(0)  are  statistically  independent  a  priori  with  Normal 

(*0.  P0>  an(*  Normal  (0,  L0)  distributions  respectively.  Zero-mean  Gaussian 
white  noise  processes  Wj,  w2,  and  v  are  independent  with  covariance 

parameters  Q,  D  and  R/m  ,  respectively.  Positive  scalars  m  and  h  and  are 
considered  parameters  of  the  problem,  with  1«  m  «l/h  and  usually  further 
restrictions.  The  objective  is  to  find  a  control  law  that  minimizes  a  scalar 
performance  criterion  of  the  form 


J  = 


T 

XfSfXf  + 


(xTAx  +  uTBu  +  2xTZu)dt 


(4) 


where  <f  denotes  prior  expectation  and  tf  is  some  specified  terminal  time.  As 

usual,  a  control  law  is  defined  as  a  decision  rule  that,  for  each  t  in  [0,  tf) , 

specifies  the  current  control  u(t)  as  a  function  of  the  current  measurement 

history  {(z(s),  s)  :  0  <  s  <  t).  A,  B,  and  Sf  are  symmetric  with  A,  Sf,  P0,  L0,  Q,  and  D 
positive-semidefinite  and  B  and  R  positive-definite.  In  the  context  of 

sufficiently  small  h,  the  components  of  all  the  matrices  as  well  as  R1  and  B'1 
are  of  order  unity.  These  matrices  may  be  time-varying  but  only  with  time- 

derivatives  that  are  also  of  order  unity  in  this  sense. 

Also,  some  further  structural  properties  are  assumed  for  the  limiting 
form  of  this  problem  when  h  =  0  (in  which  case  0  is  irrelevant).  First,  it  is 

assumed  that  (F,  G)  is  stabilizable  and  (F,  H)  is  observable  at  each  time  instant. 

To  specify  the  other  properties  and  also  for  future  reference,  it  is  convenient 

to  define  V(t)  (a  normalized  covariance  matrix  of  x  for  this  limiting  case)  and  a 
corresponding  transition  matrix  T  by 


V  =  FV  +  VFT  +  m  (Q  -  VHTR'1HV);  V(0)  =  mP0 


(5) 


and 


— =  [F  -  mVHTR_1Hl  (s)  T(s,t);  T(t,t)  =1.  (6) 

ds 

Other  assumptions  arc  made:  with  m  considered  as  a  problem-formulation 

parameter,  there  exists  a  positive  x(m)  (a  "settling  time")  such  that 
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lim  x=0andlim  mx>l, 

m  — m— 


and  for  all  sufficiently  large  m  and  t  »  xt 

( 1 )  V  is  of  order  mx 

(2)  T(s,t)  is  of  order  mxe 


for  s  >  t. 


(3)  T(s,t)VHT  is  of  order  c'(s'l)/z  J 

(4)  If  U  and  X  are  symmetric  matrix  time  functions  of  order  unity  with 

X  positive-semidefinite  and  U  and  X  also  of  order  unity,  and  if  K  and 
q  arc  the  matrix  time  function  and  vector  random  process  defined  by 


-K  =  K(F-mVHTR'1H)  +  (FT-mHTR1HV)K  +  U;  K(tf)  =  0 


and 


q  =  (F-mVHTR1H)q  +  mw3;  q(0)  =  0, 


where 

w3  is  zero-mean  Gaussian  white  noise  with  covariance 
parameter  X,  then 

(a)  TT(s,t)KVHT  is  of  order  mx2  for  s  >  t 

and 

(b)  HVKq  is  of  order  m2x5/2 

(except  perhaps  for  a  negligibly  improbable  set  of  realizations  of  q). 

Because  T(s,t)  decays  rapidly  with  increasing  (s  -  t)  >  0,  the  variations  of  V(s), 
H(s),  and  K(s)  are  negligible  until  T  is  essentially  zero;  thus,  (3)  and  (4a)  will 
also  hold  when  V,  H,  and  K  are  evaluated  at  s  instead  of  t.  Also,  since  V/m  is  the 
error  covariance  matrix  for  the  Kalman-filter  estimate  of  x,  these  assumptions 
imply  that  the  entire  vector  x  can  be  tracked  arbitrarily  closely  if  the 

measurement  noise  can  be  made  small  enough  (by  (1),  (2),  and  x-*0  as  m-»°°). 

Finding  an  optimal  control  law  exactly  is  very  difficult  for  this  case,  and 

we  seek  only  an  approximation  that  is  asymptotically  accurate  to  order 

9  9  S  /  9  0  '  1/9 

h  m  x  for  small  h  and  h  m  x  «  1  (except  perhaps  for  a  set  of 

measurement  histories  of  negligible  probability).  However,  since  x  can  be  as 
large  as  m'1/n  when  x  is  n-dimensional  (see  next  section),  control  terms  of  this 
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2 

magnitude  can  be  large  compared  to  h  under  these  conditions  if  n  >  2,  which 

is  necessary  for  their  significance  in  the  type  of  applications  described  in  the 

last  section.  The  treatment  of  the  problem  here  is  further  limited  to  finding 
the  control  law  corresponding  to  a  cost-to-go  function  that  has  the  formal 

appearance  of  satisfying  the  Bellman  equation  for  Equations  1  through  4  to 

order  h  m  x  .  This  control  law  would  be  the  desired  approximation  if  the 
equations  involved  in  the  analysis  are  well-posed — and  the  formally  negligible 
terms  are  indeed  so  in  some  appropriate  sense — but  a  rigorous  verification  of 
these  conditions  is  beyond  the  scope  of  this  report.  In  this  sense,  the  result 

obtained  below  is  only  a  plausible  candidate  for  the  control-law  approximation 
being  sought. 


HEURISTIC  MOTIVATION  FOR  RESTRICTIONS 

Although  not  readily  apparent,  the  extra  restrictions  (1)  through  (4)  in 
the  problem  formulation  are  actually  motivated  by  the  prototypical  case  of  an 

n-dimensional  system  with  x;  =  xi+1,  i  =  I,...,  n  -  1,  xn  =  order-unity  noise  w,  and 

scalar  measurement  of  xj  only.  One  way  of  estimating  the  current  value  of 

this  state  vector  for  m  »  1  (i.e.,  small  measurement  noise)  is  to  partition  a  short 
preceding  time  interval,  say  of  length  x,  into  n  equal  subintervals,  estimate 
the  value  of  xt  at  the  midpoint  of  each  subinterval  as  the  average  value  of  z 

observed  during  it,  and  then  estimate  xi+1(=  x(l^),  i  =  1,  ....  n  -  1  by  taking  an  ith- 
order  difference  of  these  Xj -estimates  and  dividing  by  the  appropriate 

multiple  of  x1.  By  standard  results,  the  error  variance  in  the  xj  estimates  in  the 
absence  of  the  process  noise  w  is  m  x"  and  that  of  the  Xj-estimate  is  of  order 
m '  x  The  error  variance  in  the  estimate  of  x;  would  therefore  decrease 

with  increasing  x  until  the  variance  of  the  change  in  Xj,  produced  by  the 
process  noise  over  that  interval  (which  increases  with  x)  is  of  the  same  order 
of  magnitude.  For  each  i,  this  happens  when  x  is  of  order  m '^n,  because  X;  is 

the  (n  -  i  +  l)st  integral  of  w,  whose  variance  over  an  interval  of  length  x  can 
be  obtained  by  standard  methods  and  is  of  order  T2(n  i)  +  1.  For  any  given  n, 
m'1/n-»0  as  m-»°o,  so  an  optimal  x  for  this  kind  of  estimation  would  indeed  be 
short  and  have  this  order  of  magnitude  for  sufficiently  large  m  (i.e.,  small 

enough  measurement  noise).  The  corresponding  error  variance  in  the 
estimate  of  Xj  would  then  be  of  order  xn+1'2i/m.  Since  the  (n  x  n)  covariance 
matrix  P  for  these  errors  must  be  positive-semidefinite,  component  Py  would 

be  of  order  xn+1','Vm. 

The  basic  premise  here  is  that  the  Kalman  filter  for  such  a  system  should 
behave  in  the  same  way  to  an  order-of-magnitude  approximation  (this  has 
been  verified  foi  some  special  cases).  In  the  notation  of  the  preceding  section 
for  such  a  system,  the  properties  assumed  there  can  then  be  established  by 
concluding  in  turn  that: 
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(1)  Vjjis  of  order  x 1  'l  J/m;  because  V  =  mP 

(2)  (VHT)j  is  of  order  x‘‘/m;  because  H  =  [1,  0,  ....  0]  in  this  case 

!of  order  x3'1  for  positive  t  -  s  of  order  x 
negligible  for  x  -  s  »  x 

in  order  for  the  transition-matrix  solution  (e.g..  Reference  6)  of  the  Kalman 
filter's  differential  equation  for  the  conditional  mean  to  operate  (as  a 
convolution)  on  the  measurement  history  in  the  manner  described  above. 

(4)  Kjj  is  of  order  m2x1+1+^;  from  the  transition-matrix  solution  of  K, 

(3),  and  the  fact  that  m  is  of  order  x‘n  in  this  case. 

(5)  var(qi)  is  of  order  m2x3'2\  from  (3)  in  the  transition-matrix  solution 

of  the  usual  differential  equation  for  the  covariance  matrix  of  q, 

and  the  fact  that  x  «  1  and  m  is  of  order  x'n. 

(6)  q{  is  of  order  mx3/2_1;  by  (5),  the  Chebychev  inequality,  and  the  fact 
that  q  is  zero-mean. 

(7)  V  is  of  order  mt;  from  (1)  and  the  fact  that  x  «  1  and  m  is  of  order  x'n. 

(8)  T(t,s)  is  of  order  mxe^'1^  for  t  >  s;  by  (3)  and  the  fact  that  x  «  1  and 

m  is  of  order  x'n. 

(9)  T(t,s)VHT  is  of  order  e^s'^/x;  from  (2),  (3),  and  the  fact  that  x  «  1  and 
m  is  of  order  x~n. 

(10)  TT(l,s)KVHT  is  of  order  mx2;  from  (2)  through  (4)  and  the  fact  that  x  «  1 
and  m  is  of  order  x'n. 

(11)  HVKq  is  of  order  m2x5/2;  from  (2),  (4),  (6),  and  the  fact  that  m  is  of 
order  x"n. 

In  this  case,  t  would  be  a  meaningful  "settling  time”  for  this  Kalman  filter. 
Although  x  can  be  tracked  arbitrarily  well  with  sufficiently  accurate 
measurements  (of  xt  alone  in  this  instance),  an  interesting  point  in  this 

connection  is  that  if  xk  were  also  unit-intensity  process  noise  for  some  k  <  n, 

then  the  Kalman  filter's  effective  "data  window"  for  estimating  xk+1 .  xn,  and 

the  variance  of  the  estimation  errors  for  these  state  components,  would 
remain  of  order  unity  as  m— »«>.  The  limiting  factor  in  estimating  these 
components  would  then  be  the  ratio  of  the  process  noise  intensities  to  each 
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other,  and  the  entire  state  vector  could  no  longer  be  tracked  arbitrarily  well 
just  by  reducing  the  measurement  noise  even  though  the  system  would  still  be 
observable. 

Finally,  it  is  reasonable  to  hope  that  the  properties  assumed  ;n  the 

problem  formulation  would  also  be  characteristic  of  a  considerably  wider  class 
of  linear  observation  systems  whose  state  can  be  tracked  arbitrarily  well  with 
small  enough  measurement  noise.  The  reasoning  is  that  the  tracking 
accuracy  in  such  cases  would  often  be  limited  by  some  similarly  configured 

subset  of  the  state  components,  one  of  which  is  driven  by  process  noise  and  is 

a  higher  derivative  of  another  that  is  measured  directly. 


STATE  ESTIMATION 

If  the  conditional  mean  and  covariance  matrix  for  the  current  composite 
state  (x,0),  given  the  currently  available  measurements,  are  partitioned  in  the 
obvious  way  as 


then  it  follows  from  standard  Kalman  filtering  results  (Reference  7)  that 


x  =  Fx  +  Gu  +  m2[P  +  h(rTET)  u]HTR'1(z  -  z);  x(0)  =  x0  ( 7 ) 

8  =  m2[ET  +  h(rTL)  u]HTR'I(z  -  z);  0(0)  =  0O  ( 8 ) 

P  =  FP  +  PFt  +  Q  -  m2[P  +  h(rTET)u]HTR  IH[P  +  h(Er)u];P(0)  =  P0  (9) 
E  =  FE  -  m2[P  +  h(rTET)'ujHTR‘1H[E  +  h(Lr)'u];  E(0)  =  0  (10) 

L  =  h2D  -  m2[ET  +  h(rTL  )u]HTR1H[E  +  h(Lr)’ul;  L(0)  =  L0,  (11) 

where 

z  =  H[x  +  h  tr(r  0uT)].  (12) 

It  will  be  convenient  to  define  the  quantities 


and 


M  =  —IE  +  h(LDu] 
h 


(13) 
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N  =  (mP  -  V)/h2  +  nmTr"LrTu,  (14) 

where  V(t)  is  as  defined  by  Equation  5,  and  the  normalized  Kalman-filier 
innovation  variable 

s  =  m(z  -  (15) 

Differentiating  Equations  13  and  14  and  substituting  from  Equations  7  through 
12  then  show  eventually  that 

x  =Fk  +  Gu  +  (V  +  h2N  +  h2(rTMT)'u|HTR^,  (16) 

M  =  (F  -  mVHTR1H)M  +  m(Lr)'u+  m[L(r-  TFT) 

+  h2(D  -  MTHTR1HM)r]u  -  mh2[N  +  (rTMT)'u]HTR'1HM,  (17) 

N  =  (F  -  mVHTR  1  H)N  +  N(FT  -  mHTR4HV)  -  m(VHTR1HrTMT 

+  MrHTR1HV)u  -  mh2[N  +  (rTMT)u]HTR'1H[N  +  (Mr)'u] 

+  miF(uTr  LrTu)  +  (uTr  LrTU)FT  -  — -(uTr  LrTu)].  (18) 

dt 


and 


L  =  h2(D-MTHTR'IHM).  (19) 

L  clearly  remains  positivc-scmidefinite  and  of  order  unity  over  time  intervals 
that  arc  short  compared  to  1/h  .  Also,  %  can  be  treated  as  a  zero-mean  Gaussian 
white  noise  process  with  covariance  parameter  R(t)  in  determining  the 
statistical  behavior  of  the  above  Kalman-filtcr  variables  (Reference  4). 

Controls  with  Wiener-process  components  will  be  used  below;  they  arc 
not  actually  differentiable  and  so  should  really  be  interpreted  as  smooth 
approximations  of  the  Wiener  process  used  formally.  Such  approximations 
would  be  unavoidable  in  practice  anyway  and  presumably  could  be  made  in  a 
way  that  amounts  o  local  averaging  over  a  short  enough  time  interval  that 
the  difference  wf  ^id  be  negligible  to  the  order  of  accuracy  retained  in  the 
analysis. 


APPROXIMATE  ESTIMATOR  BEHAVIOR  IN  A  RESTRICTED  CONTEXT 

If  h  =  0,  it  is  a  standard  result  that  the  optimal  control  law  is 

u  =  -  W  x ,  (20) 
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where  W  =  B_1(ZT  +  GTS)  with  S(t)  the  solution  of 

.  -S  =  SF  +  FTS  +  A  -  (SG  +  Z)B'1(Z7  +  GTS);  S(tf)  =  Sf.  (21) 

Since  we  are  only  concerned  with  small  h  here,  we  consider  control  laws  of 

‘  the  form 

u  =  -Wi+q,  (22) 

0  (i  1  1  /7 

under  the  further  restriction  that  h  m  x  «  1,  such  that  q  is  small  compared 
to  unity  in  that  context.  With  respect  to  sufficiently  small  h,  it  is  clear  that  the 
assumed  properties  of  the  matrices  defining  the  problem,  including  the 

stabilizability  of  (F,  G)  and  the  observability  of  (F,  H),  imply  that  W  and  W  can 

be  considered  of  order  unity  and  that  x  and  x  are  also  of  order  unity  by  this 

standard  (except  perhaps  for  a  set  of  realizations  of  negligible  probability) 

when  the  control  law  of  Equation  20  is  used  and  h  =  0. 

Differentiating  Equation  22  and  using  Equation  16  shows  that 

u  =  -(W  +  WF  +  WGW)x  -  WGq  -  W[V  +  h2N  +  h2(rTMT)' 

(q  -  Wx)]HTR’1^  (23) 

for  such  a  control  law  if  the  perturbation  control  q  is  small  and  smooth 

enough  that  q  is  negligible  in  the  following  (which  it  will  be  for  the  optimal 

control  law  approximation  derived  below).  We  now  assume  this  to  be  the  case 
and  further  that  M  and  N  are  of  respective  orders  m  x  and  m  x  with  the 
autocovariance  functions  of  the  M-components  being  negligible  for  time 
differences  that  are  large  compared  to  x.  Then  it  follows  from  Equation  19  that 
order-unity  changes  in  the  L-components  behave  basically  as  sums  of 

l/(h2m4x4)  independent  random  increments,  each  with  a  variance  of  order 

h4m8x8.  Hence  the  variance  of  such  a  sum  is  of  order  h2m4x4,  which  is  small 

compared  to  unity  under  the  restrictions  of  interest  here  because 

1  »  h2m6x11/2  =  (h2m4x4)(mx)3/2V"m.  This  means  by  the  Chebychev  inequality 

that  the  difference  between  L  and  its  prior  expected  value  L(t)  is  small 
compared  to  unity.  In  consequence,  the  substitution  from  Equations  22  and  23 

,  for  u  and  u  in  Equations  16  through  18  gives 

x  =  (F  -  GW)x  +  VHTR  (24) 

M  =  (F  -  mVHTR-'H)M  -  m(Lr)  WVHTR',t  (25) 


1  1 
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N  =  (F  -  mVHTR_1H)N  +  N(FT  -  mTR'1HV)  + 

m(VHTR'1HrTMT  +  MrHTR'!HV)W£,  (26) 

except  for  higher-order  terms,  when  M  and  N  are  of  the  assumed  orders  of 
2  6  11 12 

magnitude  and  h  m  t  «  1. 

Since  R  and  R  1  are  both  of  order  unity  and  since  Equation  24  is  the 

equation  for  x  when  h  =  0,  the  components  of  the  product  VH  must  be  of  order 
unity  because  x  is  of  order  unity  then. 

Let  q  denote  an  arbitrary  column  of  M  and  k  its  index.  Then  from 
Equation  25, 

q  =  (F  -  mVHTR'1H)q  -  mUWVHTR'^;  q(0)  =  0,  (27) 

where  U  is  the  matrix  such  that  Uy  =  CLF  )kj j-  Let  C(t)  denote  the  prior 
covariance  matrix  of  q(t).  It  follows  from  Equation  27  and  standard  results  that 

C=  (F  -  mVHTR"1H)C  +  C(FT  -  m^R^HV)  +  m2UWVHTR  1HVWTUT; 
C(0)  =  0 


and  that 


t 

C(t)  =  m2 1  T(t,s)[UWVHTR1HVWTU](s)TT(t,s)ds. 

Jo 

T 

From  the  properties  assumed  for  the  transition  matrix  T  and  the  fact  that  VH 
is  of  order  unity,  C  is  of  order  m4r3.  By  the  Chebychev  inequality,  therefore, 
the  components  of  any  column  of  M  and,  hence,  any  component  of  M  itself 
would  be  of  order  m2,c3/2  in  this  approximation,  except  perhaps  for  a  set  of 
realizations  of  negligible  probability. 


Furthermore  if  C(t2,tj)  is  used  to  denote  the  prior  expectation  of 
q(tj)qT(t2)  for  t2  >  tj  ,  it  is  a  standard  result  that 


rOt2,t,) 

dt2 


[F  -  mVHTR"1Hl(t2)  C(l2,t,);  Ca,.tt)  =  C(t,). 


and,  hence,  that  C(t2,t j )  =  T(t2,t,)C(t,).  Since  T(t2,t,)  decays  exponentially  as 
(t2  -  tj  )/t  by  assumption,  C(t2,t  j )  is  indeed  negligible  for  t2  -  tj  »  t  . 
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Now  define  X(t)  =  T(s,t)M(t)  for  t  s  s,  so 

X  =  +  T(s,t)M(t). 

3t 

From  standard  results  for  transition  matrix  derivatives  and  the  above 
approximation  for  M, 

-X  =  mT(s,t)(Lr)'wVHTR'14;  X(s)  =  M(s) 

except  for  higher-order  terms.  Integration  in  reverse  time  from  t  =  s  gives 
X(t)  -  M(t)  as  a  zero-mean  random  variable  with  component  variances  of  order 

nt4x3  because  of  the  order  of  magnitude  assumed  for  T.  Since  M(t)  itself  is  of 
order  m2x3/2,  so  is  T(s,t)M(t)  for  t  <  s  by  the  Chebychev  inequality. 

With  routine  matrix  manipulation,  the  solution  to  the  approximate 
equation  for  N  can  be  expressed  as 

t 

N(t)  =  m  I  {T(t,s)[QT(s)MT(s)  +  M(s)n(s)]TT(i,s)}'W(s)x(s)ds, 

Jo 

T  1 

where  £2  denotes  TH  R'  HV.  Because  of  the  preceding  order-of-magnitude 
results,  this  shows  that  N  is  of  order  m3x5/2. 

Under  the  restrictions  specified  above,  the  Kalman  filter  equations 

(expressed  in  terms  of  V,  M,  N,  L,  x,  and  8)  are  therefore  formally  consistent 
with  the  orders  of  magnitude  and  correlation-function  behavior  assumed  for 
M  and  N,  and  verify  them  under  the  assumption  that  the  full  equation  system 
is  well  posed  in  some  approximate  sense,  and  under  the  further  assumption 

that  t|  would  be  negligible  in  the  preceding  derivations.  As  further 
consequences,  moreover,  VH  is  of  order  unity,  T(s,t)M(t)  is  of  order  m  x  for 

s  >  t,  and  the  difference  between  L(t)  and  its  prior  expected  value  L(t)  is  small 
compared  to  unity  under  these  conditions. 


CONTROL  OPTIMIZATION 

For  S ft)  and  u  as  specified  by  Equations  21  and  22,  the  problem  here 
reduces  to  that  of  finding  an  optimal  control  law  for  the  perturbation  control  ri 
to  which  we  seek  only  an  asymptotic  approximation.  An  optimal  expected  cost- 
to-go  function  can  be  defined  consistently  in  terms  of  time  and  the  conditional 
distribution  of  x  and  9  (Reference  8).  The  Principle  of  Optimality  of  dynamic 
programming  can  then  be  applied  in  the  usual  way  (Reference  9)  to  derive  a 
Bellman  equation  for  this  function,  the  solution  of  which  equation  specifies 
the  optimal  control  law  for  tj. 
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Following  this  approach,  we  consider  a  possible  cost-to-go  function  of  the 

form 


J  =  ^rkT Sx.  +  ~  trf(S  +  Y)N]  +  h2xT  tr(A"M)  + 

2  2  m 


f(o, 


(28) 


where  S(t)  is  as  defined  by  Equation  21,  and  is  hence  symmetric  and  of  order 
unity,  and  where  Y  and  A  are  functions  of  t  only,  with  Y  symmetric  and  of 
2  3  3 

order  m~T  ,  and  with  A  of  order  nu  .  Retaining  accuracy  only  to  orders 

o  5/0  ^  25/2 

h“m  x  and  h  m  x  y  (where  y  denotes  any  product  of  q  components)  and 
using  Equations  16  through  18  (with  Equations  22  and  23  substituted  for  u  and 

u)  to  evaluate  the  required  conditional  expected  increments  (and  products 

thereof)  in  x,  M  and  N  result  in  a  Bellman  equation  of  the  form 


- -  min  p(t,  x,  M,  N,  q) 

Ot  n 


J(tf,  ...)  =  y{xTSfx  +  tr[Sf<Vf  +  h2N)]/m ) 


(29) 


where  p  is  a  certain  rather  lengthy  expression,  if  h2m6t11/2  «  1  «  m6x5/2.  This 
derivation  also  uses  the  fact  that  L  can  be  replaced  by  its  prior  expected  value 

L  in  Equations  17  and  18  in  this  context.  The  minimizing  q  in  Equation  29 
occurs  when  the  q -derivative  of  the  expression  p  is  zero,  which  is  easily  shown 
to  be  when 


> 


q  =  h2B'1tr([rHTR'1HVY  -  (AG)"]M).  (30) 

Substituting  from  Equation  30  for  q  in  Equation  29  and  collecting  terms  in  like 
powers  of  the  a  priori  random  variables  x,  M,  and  N  then  show  that  J  of 
Equation  28  formally  satisfies  the  Bellman  equation  to  order  h  m  x  '  ,  for 
h2m6x11^2  «  1  «  m6x5/2  and  S(t)  the  solution  of  Equation  21,  if 

-Y  =  Y(F  -  mVHTR'1H)  +  (FT  -  mHTR  ]HV)Y  + 

(SG  +  Z)B  l(ZT  +  GTS);  Y(tf)  -  0  (31) 

-a  =  A[F  -  GB  '(ZT  +  GTS)]  +  (FT  -  mHTR  'HV)A  + 

YVHTR1HrB1(ZT  +  GTS);  A(tf)  =  0  (32) 

f(t)  =  |  g(s)ds. 
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if  the  resulting  values  of  Y  and  A  have  the  assumed  orders  of  magnitude  and  if 
the  corresponding  perturbation  control  of  Equation  30  is  such  that  the 

contributions  of  T|  would  be  insignificant  in  the  computation  of  N  and  M  . 

Now  by  standard  properties  of  transition  matrices, 
c'f 

Y(t)  =  J  TT(a,  t)lSG  +  Z)B‘'(Zt  +  GTS)](a)  T(a,  t)  da. 

Since  T(a,  t)  decays  to  zero  by  assumption  over  an  a-interval  of  order  t  and  has 
maximum  component  magnitudes  of  order  mt,  this  integration  is  effectively 

2  3 

over  an  interval  of  order  x  and  Y  is  therefore  of  order  m  x  ,  as  assumed.  Also, 
for  the  transition  matrix  T2  defined  by 

3T2(a,t)  ,  T 

— =■ - -  [F  -  (SG  +  Z)Es  Gua)  T2(a,  t);  T2(t,  t)  =  1, 

3a  z  z 

it  is  straightforward  to  verify  that 

A(t)  =  TT(a,  t)[YVHTR'1Hr’B‘1(ZT  +  GTS)](a)  T2(a,  t)  da. 

Since  T2(a,  t)  is  order  unity  and  this  integration  is  also  effectively  over  an 
interval  of  order  x,  A  is  of  order  mx  ,  as  assumed,  by  the  unity  order  of 
magnitude  postulated  for  matrix  products  of  which  T  YVH  is  an  example. 
Finally,  the  time  derivative  of  the  derived  optimal  perturbation  control  will, 

7  7  0  A  *  <1  7  1  ' 

as  a  result,  be  dominated  by  terms  of  orders  h  m  x  M,  h  m  x  M,  h  m  x  M,  and 
h'm  x  M.  Such  terms,  as  was  assumed,  would  make  no  significant  contribution 
as  a  part  of  u  in  the  differential  Equations  17  and  18  for  M  and  N  if  h2m6x11/2  «  1 

To  summarize,  it  has  been  shown  at  least  formally  that  the  optimal  control 
law  is 


u  =  -B^fZ1  +  GTS)x  +  ti 

to  order  h2m2x5/2  if  h2m6xn/2  «  1  «  m6x5^2,  where  S(t)  is  given  by  Equation  21, 
t|  is  given  by  Equation  30  with  V(t),  Y(t),  and  Aft)  as  given  by  Equations  5,  31, 

and  32,  and  where  x  and  M  are  generated  from  the  incoming  measurements  z 
by  the  Kalman  filter  of  Equations  7  through  12  and 

M  =  (F  -  mVHTR''H)M  -  m2(Lr)  B1(ZT  +  GTS)VHTR  !(z  -  z); 

M(0)  =  0.  (33) 
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The  last  equation  is  a  sufficiently  accurate  truncation  of  Equation  17  with 
substitution  from  Equations  22  and  23  for  u  and  u. 


APPLICATIONS  IN  A  MORE  GENERAL  CONTEXT 

The  preceding  results  also  have  applications  to  a  wider  class  of  optimal 
control  problems,  with  multivariate  dynamics  of  the  general  form 

1x  =  f(x.  0,  u,  t,  v)  (motion-state) 

0  =  0  (parameters) 

|  x ( t0 ).  0  Normal  and  independent  a  priori  with  <f  (0)  =  0, 

state  measurements  of  the  form 


{/.  =  g(x,  0,  u,  t,  v) 

w,  V  independent  white  noise  processes, 
and  performance  criterion  (to  be  minimized)  of  the  form 

J  =  <f  [y(xf,  0)  +  J  X(x,  0,  u,  t)dt], 

where  t0  and  tf  are  specified  a  priori.  Corresponding  to  such  a  problem,  there 
is  the  nominal  control  problem  (easier  to  solve  because  it  is  deterministic)  of 

finding  u(t)  to  minimize 


J=  v(xf. 


X[x(t),  0,  u(t),  t]dt 


for 


x  =  f(x,  0,  u,  t,  0);  x(t0)  =  prior  mean  of  x(t0). 
Also,  a  nominal  measurement  history  z(t)  can  be  defined  as 
z(t)  =  g[x(t),  0,  u(t),  t,  0] 
when  the  minimizing  u(t)  is  used. 
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Once  x(t),  u(t),  and  z(t)  have  been  found  for  the  solution  of  the  nominal 

control  problem,  the  original  problem  can  be  solved  by  adding  to  u(t)  the 

perturbation  control  law  for  u-u  (or  5u),  which  minimizes  J-J  (or  5J).  Taylor- 
scrics  expansions  of  f,  g,  y,  and  X  can  typically  be  used  to  express  the 

dynamics,  measurements,  and  criterion  of  this  optimal  perturbation  control 
problem  in  terms  of  the  perturbation  variables  5x,  8u,  6z,  0,  w,  and  v  (and,  of 

course,  the  predetermined  time  functions  x,  u,  and  z).  Treating  the 

perturbation  variables  as  small  and  truncating  the  Taylor-series  expansions  at 

the  lowest  nontrivial  order  of  accuracy  lead  to  a  problem  of  the  familiar 
linear-quadratic-Gaussian  form  for  the  optimal  perturbation  control  law.  The 

solution  to  this  problem  is  well  known  and  can  be  added  to  u  as  a  first 
approximation  of  the  optimal  control  law  for  the  original  problem. 

If  the  Taylor-series  approximations  are  carried  out  to  one  higher  order  of 
asymptotic  accuracy  in  this  perturbation  control  problem,  the  effect  is 

typically  to  introduce  quadratic  terms  in  the  dynamics  and  measurements,  and 
cubic  terms  in  the  criterion.  Such  a  problem  can  often  be  rescaled  so  that  the 
perturbation  variables  are  of  order  unity  and  the  coefficients  of  the  added 
higher-degree  terms  become  the  quantities  that  are  relatively  small,  say  of 
order  h.  For  cases  formulated  so  that  the  parameter  vector  0  docs  not  enter  the 
measurement  equation  in  the  linear-quadratic-Gaussian  approximation, 
possibly  by  adjoining  the  observable  subspace  of  the  original  "parameters"  to 
x,  and  with  dim(x)  >  dim(z),  this  higher-order  problem  typically  can  be 
transformed  further  to  one  with  dynamics,  measurements,  and  criterion  of  the 
fo  rm 

y  =Fy  +  Gu  +  tr(AyyT)  +  2tr(8  yuT)  +  GWN(Q),  (34) 

z  =  Hy  +  tr(<byyT)  +  2Htr(flyuT)  +  GWN(R),  (35) 

and 

J=  ^{7  yf  Syf+  j  yj  tKnfyfyj)+  f  tcTy  +  yTZu 
I  3  Jt0 

+  y(yTAy  +  uTBu)  +  yT  tri'Fuu1)  +  uTtr(YyyT) 

+  — yT  tr(IyyT)  +  -^-uT  tr(ZuuT)]dt}, 

3  3  (36) 

where  c  and  the  three-way  matrices  are  of  order  h  and  B'1  and  all  other 
quantities  arc  of  order  unity,  where  u,  z,  and  J  are  now  used  in  place  of  5u,  8z, 
and  8J,  where  "GWN"  denotes  Gaussian  white  noise  with  the  indicated 
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covariance  parameter,  where  y  denotes  the  composite  vector 
F.G,  H,Q,S,Z,  and  A  have  the  forms 

H  =  [H  :0],  Q  = 


8x 

0 


and  where 


F  = 

F  0‘ 

,  G  = 

G 

.0  0. 

.0. 

S  = 

sf  o' 

,z  = 

z 

.0  0. 

.0. 

Q0 
0  Oj 


,  and  A  = 


A  0 
0  0 


for  this  partitioning  of  y.  If  R  1  is  also  of  order  unity,  the  optimal  control  law 
for  this  problem  can  be  approximated  up  to  terms  which  are  asymptotically  of 

order  h  (sec  appendix  for  details),  and  this  approximations  retains  this 

accuracy  as  R->  0  if  ft  =  0  (Reference  10).  Furthermore,  the  Taylor-series 

truncations  in  the  problem  formulation  don't  justify  order-h  accuracy  in  the 
solution  anyway. 

For  nonzero  ft,  however,  the  analysis  of  the  preceding  sections  shows 

that  if  R  =  R/m  ,  and  R,  F,  and  H  have  the  properties  assumed  there  with 
h2m6Tll/2  «  1,  then  the  optimal  control  law  can  also  have  another  component 
(Equation  30)  of  order  h2m2x5f2.  This  latter  component  can  be  large  compared 

to  h  under  these  conditions  and  is  therefore  significant  at  the  level  of 
accuracy  of  the  problem  formulation.  In  this  case,  this  extra  control 
component  can  simply  be  added  to  the  (order-h)  control  law  approximation  of 

the  appendix  to  approximate  the  optimal  control  law  to  order  h  m  z  . 
Absorbing  h  and  m  into  the  definitions  of  T  and  R  in  earlier  sections  allows 
this  added  control  component  to  be  expressed  as 

T1  =  B  itr{[rHTR'1HPY  -  (AG)"]  M),  (37) 


where  f  is  the  p  x  m  x  n  three-way  matrix  such  that 


n.j.k  ~  ^j,k,i  +  n 


n  =  dim(x) 

m  =  dim(u) 
p  =  dim(0)  , 


where  P(t),  Y(t),  and  A  arc  the  solutions  of 
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-Y  =  Y(F  -  PH^H)  +  (FT+HTR4HP)Y  +  (SG  +  Z)B~'(ZT  +  GTS); 

Y(tf)  =  0  (39) 

-A  =  A[F  -  GB_1(ZT  +  GTS)]  +  (FT+HTR'IHP)  A  + 

YPHTR',Hr  B''(ZT  +  GTS);  A(tf)  =  0,  (40) 

with  S(t)  as  defined  by  Equation  21,  and  where  M  is  generated  by 

M  =(F  -  PHTR  !H)M  -  (LOB’^Z7  +  GTS)PHTR1(z  -  Hi); 

M(t0)  =  0  (41) 

with  L  denoting  the  conditional  covariance  of  0  (as  computed  in  any  of  the 
preceding  approximations  of  the  optimal  perturbation  control  law). 

This  added  control  component  has  the  basic  form  of  an  irregular  dither 

term  that  is  produced  as  a  linear  combination  of  the  outputs  (M)  of  a  high- 
frequency  linear  system  (Equation  41)  whose  input  is  the  Kalman-filter 

innovation  variable  (z  -  Hx).  In  specific  cases,  moreover,  this  dither  can 
become  an  even  more  significant  component  of  the  optimal  law 
approximation,  the  result  of  particular  order-of-magnitude  properties  beyond 
what  was  assumed  in  establishing  the  more  general  results  above.  A  common 
example  is  that  of  regulating  the  position  of  a  scalar  undamped  second-order 
system  when  only  the  position  is  measured  and  only  the  velocity  is  affected  by 
the  plant  noise.  This  means  that 


F  = 

0  1 

,  G  = 

0 

,  H  =  [1 :0],  Q  = 

0  0 

,  A  = 

1  0 

-0  0_ 

.1. 

LO  qj 

.0  0. 

in  the  notation  above.  In  that  case,  it  can  be  shown  that  x  =  1/Vm"  and  that 
optimal  control  law  approximation  of  the  preceding  paragraph  is  actually 
accurate  except  for  terms  of  order  h  as  long  as  mhx  «  1.  The  dither  term  is 

still  of  order  h2m2t5/2  in  this  case  so  it  will  even  be  large  compared  to  h  under 

these  conditions  if  m'1/2  »  h  »  m'3/4.  Thus,  the  initial  parameter  uncertainty 
can  be  large  enough  in  such  a  case  without  being  too  large  to  destroy  the 
formal  validity  of  the  result  that  the  dither  control  is  the  next  level  of 
asymptotic  approximation  to  the  optimal  control  law  after  the  "linear- 
quad  rati  c-Gaussi  an"  approximation. 

Figure  1  shows  the  numerical  significance  of  adding  such  a  dither 
control  to  an  adaptive  missile  autopilot  (attitude  control  law),  which  is  an 
example  of  this  particular  type  of  control  problem  with  two  parameters 
altogether.  These  parameters  approximate  the  effects  of  missile  speed  and 
altitude.  The  details  of  this  example  arc  described  in  Reference  11  with  the 

unfortunate  exception  of  an  error  in  computing  the  dither  control  according 
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to  Equations  37  through  41  and  the  approximation  of  <f  [exp(0j)]  by  exp|<?  (0;)] 

instead  of  exp[<f  (0t)  +  yvar(0j)].  Briefly,  it  is  a  case  of  planar  motion,  and  the 

objective  is  to  make  a  missile's  acceleration  lateral  to  its  flight  path  match  an 
arbitrary  "commanded"  value  (which  typically  varies  with  time).  This  lateral 
acceleration  is  produced  primarily  by  the  body  lift  that  results  from  orienting 

the  missile  at  an  angle  to  its  flight  path,  and  this  orientation  is  controlled  in 
turn  by  tail  fin  deflections.  The  controller  measures  actual  and  commanded 
lateral  acceleration,  and  selects  the  current  tail  fin  deflection.  The 
optimization  criterion  used  for  the  control  law  is  one  for  which  the  optimal 
control  makes  the  acceleration  track  the  commanded  value  as  a  critically 
damped  second-order  system  under  idealized  conditions  (missile  dynamics  obey 
a  simplified  model,  and  parameter  values  arc  known  precisely).  Figure  1 
shows  the  result  of  using  the  indicated  control  laws  in  a  more  realistic 
computer  simulation  of  the  missile  and  aerodynamics.  The  commanded 
acceleration  in  each  case  here  is  the  indicated  scries  of  1-scc  steps.  This  is 

followed  closely  under  the  nominal  flight  condition  (Mach  2  at  20,000  ft),  both 
with  and  without  the  dither  control.  At  Mach  1.5  and  60,000  ft.,  however, 
adding  the  dither  control  caused  the  autopilot  to  adapt  more  quickly  to  this 

nonnominal  flight  condition  (i.e.,  to  nonzero  values  of  the  parameters).  The 
parameter  changes  corresponding  to  this  difference  in  flight  condition  are 

actually  far  beyond  what  has  been  justified  theoretically  for  using  this  dither 
to  approximate  an  optimal  control  law  more  accurately.  Instead,  they 
correspond  to  prior  parameter  variances  that,  to  test  the  robustness  of  this 
approximation,  were  made  as  large  as  would  allow  the  autopilot  response  to 
remain  empirically  reasonable  at  the  nominal  flight  conditions. 


acceleration  (g's) 
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without  probing  control  1  actual  accelaration  at 

with  probing  control  J  Mach  1 .5,  60,000  ft. 

commanded  acceleration 


FIGURE  1.  Autopilot  Response. 
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Appendix 

CONTROL  LAW  DETAILS 

The  order-h  approximation  of  the  optimal  control  law  for  Equations  34 
through  36  was  developed  for  the  most  part  in  Reference  10.  Including  the 
remaining  terms  in  that  development  shows  that  this  approximation  in  the 
notation  of  Equations  34  through  36  is 

u  =  -B'1(Zr  +  GrS)y  +  B'1  tr{[2(SG+  Z)B  V  -  T  -  0S 
-S0T  -  (SGt-Z)B'1HB'1(ZT  +  GTS)  +  (nG)"]y  y  T} 

-B"1  {c  +  GT4>  +  tr[TP  +  (S  +  Y)(P0  +  0TP)  +  PHTR‘,HPY}> 

where  y  and  P  arc  generated  from  the  incoming  measurements  z  by 


y  =  (F  +  0  u)y  +  (G  +  0Ty)u  +  tr[A (y  yT  +  P)] 


+  P(H  +  2Hfi  u  +  20 "y)TR_1(z  -  z);  y(t0)  = 


£  [x(t0)l 
0 


and 


P  =  (F  +  2A  y  +20  u)P  +  P(F  +  2A  y  +  2©  u)T  +Q 

-P(H+  2Ho"u  +  20 "y  )TR'1(H  +  2Ho"u  +  2<b”y)P 


+2[AHT  +  (POP)’]R'1(z  -  z);  P(t0)  = 


var[x(t0)l  0 
0  var(0)J 


and  where 


z  =  Hy  +  tr[<t>  (y  y  T  +  P)] 

and  A  (t),  S(t),  Y(t),  ri  (t),  and  <J>(l)  are  the  solutions  of  the  differential  equation 
system 


V  =  FV  +  VFT  +Q  -  VHTR' 1  HV;  V(t0)  =  P('0) 
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a  =  r  +  r'  +  r";  r  =  (F  -  vhtr‘  1  h>a  +  vav  -  (Vov)R' 1  hv;  A(t0)  =  o 

-S=SF  +  FTS  +  A  -  (SG  +  Z)B  1(ZT  +  GTS);  S(tf)  =  S 

-Y  =  Y(F  -  VHTR'rH)  +  (FT  -  HTR  1  ~HV)Y  +  (SG+ Z)B_1(ZT  +  GTS);  Y(tf)  =  0 

-n  =  i  -  [(SG  +  z)b1zb1(zt  +  gtS)]'b'1(zt  +  gtS) 

+r  +  r'  +  r";  r  =  n[F  -  gb  _1(ZT  +  gtS)]  +  a's 
-s©  b1(zt  +  gtS)  -  (SG  +  z)B~’e'Ts  -  r  at 
+(SG  +  Z)B~1(H'  +  v")B'1(ZT  +  GTS);  n(tf)  =  nf 

-4)  =  [FT  -  (SG  +  ZJB-'g'1'^  -  (SC  +  Z)B_1  (c  +  tr[TV 

+(S  +  Y)(V0  +  eTV)  -  V(QHTR‘ 1  H  +  HTR  1  HQT)VY] } 

+S  tr(VA)  +  tr[VE  +  VHTR1HVn  +  (S  +  Y)(Va’  +  a"v) 
-V(HTrV  +<DR‘1H)VY];  0(tf)  =  tr[nfV(tf)]. 
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